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Transition from stable orbit to chaotic dynamics in hybrid systems of Filippov type 
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We demonstrate on a representative example of a planar hybrid system with digital sampling a 
sudden transition from a stable limit cycle to the onset of chaotic dynamics. We show that the 
scaling law in the size of the attractor is proportional to the digital sampling time r for sufficiently 
small values of r. Numerical and analytical results are given. The scaling law changes to a nonlinear 
law for large values of the sampling time r. This phenomenon is explained by the change in the 
boundedness of the attractor. 
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The control, design and analysis of many systems rel- 
evant to real world applications involves understanding 
the interaction between continuous and discrete dynam- 
ics. For example, the automated control of a car moving 
on a road is implemented by digital computer but the mo- 
tion of a car is continuous in time [l| . Hence the design 
of the control of such a system needs to take into account 
the effects of the interaction between the continuous and 
discrete dynamics. Other examples include the control of 
the motion of digitally controlled machines, for instance 
in robotics 0,13, [jj- These types of control systems are of- 
ten referred to in the control literature as hybrid control 
systems !• In 0,i it has been shown that the digiti- 
zation of the spatial structure by the controller induces 
micro-chaotic transient dynamics. Effects of digitization 
on the stability of the solutions have been considered in 
[f| 0]> an d in [H the existence of different types of at- 
tractors in a simple model of a delta-modulated control 
system has been shown. Here we consider another as- 
pect of digitization in hybrid systems. We assume that 
the input to the controller is delivered at discrete times, 
separated by a constant r > 0, and show that for arbi- 
trarily small t the system can exhibit chaotic dynamics. 
Moreover, there are scaling laws relating the maximum 
distance of the chaotic attractor due to digital sampling 
from the simple periodic attractor of the continuously 
sampled system. It is linear for sufficiently small values 
of r but at larger r there is a change in the properties of 
the boundedness of the chaotic attractor and the scaling 
becomes nonlinear. We will illustrate this effect by con- 
sidering simple Filippov systems for which the evolution 
of a variable x in some region D C l n is determined by 
the equations 



±{t) 



Fx(x(t),fi) if H(x(t),(jt)>Q 



(1) 



f 2 (x(t),n) if H(x(t),n)<0, 
where F\, F 2 are sufficiently smooth vector functions 



and H(x(t),fj,) is some smooth scalar function depend- 
ing on the system states x £ K n , and parameter /i £ K m ; 
t 6 1 is the time variable. The boundary £ on which 
H(x, /i) = is assumed to be a hyperplane which divides 
the region D into two subspaces, G\ (for H(x, fi) > 0) 
and G 2 (for H(x, fi) < 0), in which the dynamics is 
smooth and continuous. There may be extra specifi- 
cations which determine the motion across or in E, for 
example if there is sliding motion as is the case in our 
example, but these are standard to include. 

In ([T|) the control of the switching between the two 
systems across £ is instantaneous. The modified hybrid 
Filippov systems we study are obtained by assuming that 
the control function H is evaluated at discrete times kr, 
k = 0,1,2,..., for some constant r > 0, and so the 
decision to change evolution equation can only occur at 
these discrete times. Thus for k = 0, 1, 2, 3, . . . wc define 
a discrete variable ik by 



1 if H(x(kr),n) > 
//, i = { 2 if H(x{kr),fi) < 
ik if H{x{kr), fx) = 0, 



(2) 



with io = I (arbitrarily chosen) so that ix is always well 
denned, and replace the evolution |T]) by 

x(t) = F lk (x(t),fjL) if (k - l)r < t < kr. (3) 

Note that this system excludes the possibility of sliding 
motion, that is a motion within the discontinuity set S. 

The example which we will consider in the remainder 
of this paper is planar. Set 



p _ i -OLX\ — UJX2 + xi{x\ + x\) 

1 ' ujxi — ax 2 + x 2 (x\ + x|), 



Fn 



b. 



(4) 
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H(x) = x 2 — n, and h — kr (k — 0, 1, 2,- ■ •) where 
a, uj, a, 6, fi and r are some chosen constants (system 
parameters). The vector field F\ is the normal form for 
a simple subcritical Hopf bifurcation, with a stable focus 
at the origin and an unstable periodic orbit with radius 
yfa. if a > 0. For appropriate choices of fi, a and b 
the vector field F 2 can be used to create a stable periodic 
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(a) (b) 
FIG. 1: Asymptotic trajectories in Filippov system Q with- 
out digital sampling (a), and with digital sampling r = 0.01 
(b). 

orbit for the Filippov system {l} as shown in Figure Ufa). 
The stability is derived from the fact that part of the 
cycle lies on S, and this segment of the orbit is called a 
sliding segment |9|. This stable cycle may coexist with 
the unstable cycle of the vector field Fi . 

For all the following numerical computations we set 
w = 15, a = 1. In this case the vector field unstable limit 
cycle of Fi is centered at the origin and has radius 1. 
The vector field F2 is assumed constant; we set a = — 1 
and 6=1. Finally let us set /x = 1.1 so that £ does not 
intersect the unstable limit cycle of Fi . For these param- 
eters (QJ with vector fields given by ([4]) admits the stable 
limit cycle with sliding already referred to as well as the 
unstable cycle of Fi (the main restriction in the choice 
of the parameters is that given //, the angular velocity u> 
in F\ is large enough for the sliding cycle to exist). 

Let us now increase the sampling time r from using 
(H}. As Fig.rjTb) shows, there is an apparent thickening 
of the attractor: the stable limit cycle no longer exists 
and we observe an onset of more complex asymptotic 
dynamics. It turns out that this complex dynamics is 
chaotic and so there is a transition from a stable orbit to 
a chaotic attractor due to an introduction of the sampling 
process. 

We will prove this in two parts. First we show that 
if t > is sufficiently small then there is a compact set 
which solutions cannot leave, and hence which contains 
at least one attractor. This part of the proof is based 
on showing monotonic crossing of the local transversal 
by a trajectory, similarly as in the Poincare- Bendixson 
Theorem [io| . We present this argument in fairly gen- 
eral terms below so that the extension to similar systems 
is clear. Second we show that any solution in this com- 
pact set (and hence the attractor itself) has a positive 
Lyapunov exponent. 

The compact invariant region is annular, and its inner 
boundary is the unstable cycle of -Fi with radius one. 
To begin the construction of the outer boundary, note 
that the sliding segment of the (true) Filippov system 
terminates at the point x 9 = x?,) which is where the 
solution of Fi is tangential to the surface S as shown 
in Figures [2] and [31 {x 9 2 is close to for a close to /i 2 ). 
Now, let (f>i(u,v,t) denote the flow generated by Fi, i.e. 
solutions of -Fi at time t with initial condition (u, v) and 




FIG. 2: Schematic representation of the bounding regions 
along E and the boundaries Ur and Vr in the case (a) when 
x > x, and (b) x < x. 



let </>[(«, v) — {(p,q) I {p 1 q) — <j)i(u,v,t) for some < 
t < t}. Define V to be the set of points (xi,X2) which 
can be reached from a point (u, v) S GiUE, with v > —R 
within time r, and whose trajectory intersects Gi in time 
r, i.e. the set of points (x±, X2) € G2 that are reached 
from G\ U £ within time r by following <pi . Finally, let Vr 
be the right hand boundary of V, i.e. (u,v) £ Vr such 
that if (u' , v) E V then u > vl . Along most (and in some 
examples possibly all) of its length Vr will be the time r 
image of points on £, but close to x 9 this might not be 
the case. Vr therefore represents a boundary which no 
orbit which starts in G\ above —R can cross within time 
t under the flow <f>\. The choice of R is determined by 
later considerations, it needs to be large enough to allow 
the argument to close up below - numerical experiments 
show that R = 2 is sufficiently large here. 

Now consider the effect of the flow generated by F2 
to points in G2 to the left of Vr. The trajectories are 
straight lines with slope —1, and as £2 = — 1 the furthest 
to the left that an orbit from G2 can reach in time t has 
x\ = /i — t. Let U be the union of Vr and the set of 
points on straight lines of slope —1 from (u, v) € Vr with 
«>/j-Tto/i-r. Finally let Ur be the right hand 
boundary of U. Note that Ur must be connected. 

To summarise: by construction, in time r, no solution 
of Fi above the line X2 = —R can move to the right of Vr 
(which is on the left of Ur or equal to it at places), and 
under F2 all such orbits remain to the left of Ur until 
they return to G\. 

Let x be the highest point in Ur with x± = fi and 
x the highest point in Ur with x± = /i — r. Let x be 
the first intersection of the solution of Fi through x with 

X\ — [l — T. 

If X2 > X2 then the outer boundary of the bounding 
region is the trajectory through x under Fi until it hits 
Vr for the first time (see Figure E^a)), a horizontal line 
segment from Vr to Ur, and then Ur back to x. Note 
that this requires R to be large enough so that the tra- 
jectory does hit Vr. If not then a larger R needs to be 
chosen. 

If X2 < X2 let x' be the first preimage of x on E un- 
der Fi, and note that this will lie above x. Then the 
outer boundary of the bounding region is the trajectory 
through x' under Fi until it hits Vr for the first time, 
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FIG. 3: Schematic representation of the bounding region Xb 
for small values of the sampling time r. 



a horizontal line segment from Vr to £/r, and then Ur 
back to x (see Figure Hfb)). As before R needs to be 
large enough for the connections to work. 

In either case we will have created a compact region 
which no trajectory can exit from, and hence the annular 
region contains at least one attractor. 

Consider now a sufficiently small r. Define T T F as the 
image of T under the action of 4>i for time r. Let x p S T 
be the pre-image of the point at which T F crosses T in 
the neighborhood of x 9 as shown in Fig. |^J Define Xb 
to be the set of initial points x in the neighborhood of x 9 
such that for any x € Xb a trajectory generated by Fi 
evolves through x 9 . For r sufficiently small Tip is nearly 
tangent to T and there exists a point on T F , say x b , such 
that the trajectory starting at x b crosses T at some point 
x c below x p , and the time of evolution from x b € T F to 
Xb is r (see Fig. [3j> - Therefore a trajectory starting 
at x b must switch to <j)i to the left of set Xb or on Xb- 
Moreover since x p lies above x c no trajectory generated 
solely by F\ can lie to the right of Xb - penetrate G2 
and return to G\ without switching to 2 - We further 
note that a trajectory rooted at any point within the 
region bounded by £, T F , and the line segment joining 
Xb with x c switches to the vector field F\ in a region to 
the left of Xb - the time of evolution from any point 
in this region to reach some point in G\ to the left of 
Xb is less than r. It then follows that for sufficiently 
small r, Xb is a bounding set for the attractor. To find 
the bounding set to the left of T we note that the set 
of points furthest to the left of T, which can be reached 
by a trajectory generated by F2, has co-ordinates (// — 
r, x%). This implies that along the x\ co-ordinate, in the 
neighborhood of x 9 , the difference between the largest 
and the smallest values of x\ on the attractor is r. 

For larger values of the sampling time r point x p might 
no longer lie above x c , and there exist trajectories in the 
neighborhood of x 9 solely generated by F\ that lie to 
the right of Xb- In this case the difference between the 
largest and the smallest values of x\ on the attractor is 
different from r. Therefore, we expect to see a change in 
the size of the attractor as a function of the sampling 
time, that changes from a linear law for small r, to a 
different scaling not linearly proportional to r. 

To see that the attractor born out of the stable cycle 



for non-zero r is characterized by a positive Lyapunov 
exponent we consider the determinant of the linearized 
map that maps the neighborhood of x 9 onto itself. Hence, 
we consider the determinant of the matrix composition 
of the solutions of the variational equations for the flows 
0i and <f>2- Since the vector field F2 that generates the 
flow 4*2 is a constant vector field then 



det 



dx 



where t is the time corresponding to the evolution fol- 
lowing <p2, and = <I> 2 (i) is the fundamental solution 
matrix corresponding to the flow <f>2- 

Using the explicit expression for the flow function <j>i in 
the polar co-ordinates we can compute the fundamental 
solution matrix $1, corresponding to the flow (f>i : 



dip, 9) 



I 



with j = 0, 1, 2 • • • that correspond to the times of evo- 
lution following (pi after j— th switching from the flow </> 2 
to (j>x . At j = we initialize the evolution from the neigh- 
borhood of x 9 , and pj denotes the radius from the origin 
at the j— th switching instance. Finally /j(0, pj) = 1 and 
fj are monotonically increasing functions of tj. 
Therefore, 



det $i(ij) = det 



dip, 0) 



> 1, 



and J2j tj > u ( n °te that the flow follows <\>\ from x 9 
until the first intersection with T for the amount of time 
greater than tt/lu). The resulting determinant of the 
composition of the fundamental solution matrices cor- 
responding to the flows 4>\ and (j>2 is n"~Q det <&i(tj) x f 
where n is the number of switchings that are required for 
the system trajectory to reach the neighborhood of x 9 . 
We should note here that formally we should introduce 
a co-ordinate transformation to the flow 02 and use the 
fundamental solution matrix corresponding to the vector 
field F2 in polar co-ordinates. However, since the flow 
<f>2 is a constant flow the determinant of $2 is always 1 
regardless on the co-ordinate set. 

Therefore, we are only interested in H^Tg 1 det 
Since this product is greater than f and the determinant 
of a matrix is the product of its eigenvalues we conclude 
that there is at least one eigenvalue of n™ = r o 1 <f > 1 (t J ) which 
is characterized by the magnitude greater than one. This 
eigenvalue is the exponential of the Lyapunov exponent. 

Let us determine the size of the attractor, measured 
as the distance from the largest to smallest values which 
the attractor attains on a Poincare section defined on a 
set {x2 = 0,-1.3 < x\ < —I}. For sufficiently small 
values of the sampling time r the size of the attractor 
along the x\ co-ordinate is proportional to r around the 
point x 9 . However, the size of the attractor is measured 
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FIG. 4: Size of the attractor versus the sampling time r. The 
dashed line refers to theoretical predictions of the size of the 
attractor. 



on section {x 2 = 0,-1.3 < x\ < —1}. Therefore, we 
have to determine the expansion of the attractor along 
the flow after time ir/to (which is the time required to 
map the points from the neighborhood of x 9 using flow 
4>i onto our chosen Poincare section). This expansion 
is captured by the non-trivial Floquet multiplier of the 
fundamental solution matrix. To find this multiplier we 
use the explicit solutions of the differential equations that 
define F\ in polar co-ordinates. We find that 



Pit) 



(1 — ap ) exp(2ai) — 1 



= 0q + ujt. 



Differentiating p(t) with respect to po, and after substi- 
tuting for a = 1, t = 7r/15, and for po = 1.1 we get 



dp _ 
lpo~ ~ 



1.8084. In Fig. [4] using the logarithmic scales we 

are depicting how the size of the attractor scales against 
the sampling time r. The dashed diagonal line refers to 
the linear scaling proportional to r obtained using the 
above theoretical prediction. We can see that the theo- 



retical prediction coincides with the numerical results for 
small values of r. For r sufficiently small the attractor is 
bounded by the trajectory leaving S at x 9 . 

We notice that the increase in the value of the sam- 
pling time t above r = 0.01 results in the growing dis- 
crepancy between the numerical and theoretical values. 
This comes from the fact that the attractor is no longer 
bounded by the trajectory leaving S at x 9 . Other effects 
such as resonances between the sampling time r and the 
rotation u> of the flow </>i produce other local variations 
in the scaling law visible in Fig. [4j 

In conclusion, we have shown on a planar example that 
the digital sampling applied to the decision function in 
Filippov type systems leads to the onset of chaotic dy- 
namics. Further on we have shown that for sufficiently 
small values of the sampling time r the size of the chaotic 
attractor scales linearly with the sampling time r. The 
mechanism that leads to the onset of chaos is triggered 
by the expansion of the volume of phase space produced 
by the flow </>i combined with the re-injection triggered 
by the subsequent application of the switchings along the 
manifold E. The scenario observed in our model example 
will be present in a larger class of systems with swith- 
ings. The essential ingredients of these systems will be 
the presence of a stable cycle (when no digitization is ap- 
plied) and the presence of an expansion of a volume of 
phase space in the presence of digitization. 

The onset of chaotic dynamics triggered by this mech- 
anism is similar to an abrupt transition from a stable 
periodic orbit with sliding to a small scale chaotic dy- 
namics that might occur in Filippov type systems under 
an introduction of an arbitrarily small time delay in the 
switching function 11]. On the practical side, in spite 
of the fact that these oscillations are micro-chaotic, they 
can can be highly harmful to control elements, induce ex- 
cessive wear to machine tools [■§ , and therefore it is im- 
portant to understand the mechanisms that might trigger 
this type of complex dynamics. 
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